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INTRODUCTION 

Analytical  models  for  fires  and  their  associated  environment  may 
be  useful  tools  for  predicting  and  studying  the  dynamics  of  flow  sur- 
rounding a source  of  disturbance  and  for  studying  the  production  and 
distribution  of  pollutants  that  result  from  combustion,  photochemical 
processes,  etc.  These  models  may  also  prove  to  be  useful  tools  for 
long-range  planning,  forecasting,  and  management  of  the  environment. 

It  is  emphasized  that  analytical  and/or  physical  modeling  is  a diffi- 
cult task.  It  will  not  serve  as  a substitute  for  good  experimental 
measurements,  but  it  is  a prerequisite  to  good  experimental  planning. 
Physical  models  are  readily  employed  in  systems  analysis  to  study  the 
overall  response  to  a specified  perturbation,  simultaneously  taking 
into  consideration  the  combined  effects  of  many  variables  and  parame- 
ters. The  state  of  the  art  of  fire  and  environmental  flow  dynamics  is 
not  satisfactorily  developed  at  the  present  time.  This  is  due  to  the 
lack  of  pertinent  information  that  could  be  developed  through  numerical 
simulation,  and  also  due  to  the  lack  of  sufficient  detailed  observa- 
tions in  the  field. 

THE  GENERAL  MODEL 

Flow  field  equations  are  presented  as  time-smoothed  equations  of 
rhange  for  an  incompressible  fluid,  and  subsequently  reduced  to  sim- 
plified two-dimensional  flows.  The  present  state  of  the  art  allows 
these  equations  to  be  extended  to  the  2-D  fire  front  as  well  as  the 
axisymmetric  fire  plume.  A quasi  steady-state  flow  field  results 
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when  the  time  derivative  equals  zero  and  does  not  appear  In  the  gen- 
eral equations,  but  parameters  and  boundary  conditions  may  be  time- 
d'-endent. 

The  formulation  will  be  focused  on  turbulent  flows;  however,  the 
equations  for  laminar  and  turbulent  flows  can  be  reduced  to  one  equiva- 
lent set  of  partial  differential  equations  when  the  effective  transport 
properties  are  employed  [see  reference  (2)]. 

The  macroscopic  description  of  the  state  of  a fluid  system  is  de- 
fined by  the  number  densities  of  the  various  chemical  species,  the 
mass  average  (or  stream)  velocity,  and  the  temperature.  These  variables 
are  related  to  quantities  that  are  conserved  in  molecular  collisions; 
i.e.,  the  masses  of  individual  molecules,  the  momentum,  and  the  energies 
of  colliding  molecules.  For  each  of  these  summational  invariants  there 
is  an  associated  equation  of  change;  i.e.,  the  equations  of  continuity 
of  individual  species,  the  equation  of  motion,  and  the  equation  of  energy 
balance.  Global  or  overall  continuity  is  also  preserved.  The  general 
form  of  these  equations  as  applicable  to  multicomponent  and  turbulent 
flow  problems  is  as  follows  [see  reference  (3)]: 


- (v  C.)  + (v  CA)  + |-  (v  Cjl 
L 3x  x A 3y  y A 3z  z A J 


+ 


AX 


’2  2 2 

3 C.  3 CA  3C. 

A A A 

“ 2 + . 2 .2 
3x  3y  3z 


' k“c*n 
n A 


+ 


(1-1) 


-4- 


3C, 

at 


- ■ - [h  <VV  * I7  (vyCN>  + h <\V] 


+ A 


NX 


32C 


N 


32C, 


N 


9x 


3y 


a2c. 


9z 


- k"C  n 
n N 


(1-2) 


wheis  C represents  concentration;  v ■ (v^,  v^,  v^)  is  velocity;  £ is 
the  diffusion  coefficient;  X represents  any  species  from  A to  N;  and 
k*  is  the  reaction  rate  constant.  For  simplification,  it  will  be 
assumed  that  k"  is  independent  of  position  and  is  s function  of  tem- 
perature only.  The  above  equations  are  for  the  continuity  of  species 
A,  B,  ...,  N,  which  are  assumed  to  be  generated  and/or  destroyed  by  an 
n-th  order  chemical  reaction.  In  turbulent  flows,  C^,  Cg,  ...,  will 
be  rapidly  oscillating  functions  of  time.  It  is  then  convenient  to  re- 
place these  quantities  by  the  sum  of  time-smoothed  values  C^,  Cg,  ..., 
and  turbulent  concentration  fluctuations  C^,  Cg,  ...,  C^.  Thus 


ca  - ca  + c; 


CN+  Si 


Substituting  and  v^  + v^  for  and  v^  respectively,  one 

obtains,  after  time  averaging,  the  following  forms  of  the  continuity 


equation: 
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The  underlined  expressions  represent  additional  terns  due  to  time 
smoothing  of  the  species  continuity  equation. 

The  time  smoothed  equation  for  continuity  of  the  N-th  species  can 
be  written  in  vector  notation  as 


The  symbol  J is  defined  as  J r “d  the  of  °Per*tor  1*  defined 

D 3 — 

as  — * + v • 7.  The  superscript  (t)  represents  a time-smoothed 

value,  and  (1)  represents  a departure  from  this  value. 

By  similar  considerations  the  time-smoothed  equations  of  motion, 
energy  transport,  and  global  continuity,  have  been  shown  to  be  as 
follows  [see  Bird,  et  al.,  reference  (1);  also  reference  (3)]: 
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(V  • v)  - 0 (4-1) 

In  the  above  equations,  and  are  the  Reynolds  stress  tensors; 

^ and  q^  are  equal  to  -kVT  and  -k7T  ’ respectively,  where  k is 

»v 

coefficient  of  thermal  conductivity;  q ■ q^  + q^;  and  <p  ^ and  <p 

>v  ^ r+j  V V 

represent  expressions  for  the  viscous  dissipation  function  and  the  time- 
smoothed  dissipation  function  respectively. 

The  above  set  of  equations  is  satisfactory  for  specifying  t\t 
composition  and  state  variables  in  turbulent  reactive  flows.  This 
set  of  equations  constitutes  a basic  physical  model  suitable  for  describing 
the  turbulent  flow  fields  surrounding  a source  such  as  fire  and  associated 
environment.  However,  it  is  highly  probable  that  radiation  is  an  important 
mechanism  of  heat  transfer,  and  accordingly,  a model  for  radiation  is 
desirable. 

The  heat  transfer  due  to  radiation  can  be  accounted  for  as  follows 
[reference  (4)]: 


R - A 


(5-1) 
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where  R is  the  rate  of  radiant  emission  per  unit  volume,  A is  the  rate 

of  energy  absorption  per  unit  volume,  and  r is  the  generalized  coordinate. 

For  a given  system  one  can  define  an  absorption  coefficient  ^ with  di- 
2 -1 

menslon  L M for  radiation  of  frequency  v.  If  R^(r)  is  the  energy  emitted 
per  unit  volume  per  unit  time  in  unit  frequency  range  v,  and  q^(r)  is  the 
radiation  energy  flux  in  this  frequency  range,  then 
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Equation  (t>-l)  can  be  integrated  for  two  special  cases:  opaque 

materials  (flame  fronts  with  large  amounts  of  smoke),  in  which  absorption 
or  radiation  is  so  intense  that  the  radiation  is  in  local  thermodynamic 
equilibrium;  and  transparent  materials,  in  which  there  is  little  self- 
absorption  of  radiation.  If  the  absorption  coefficient  can  be  taken  as 
constant. 
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For  weak  absorption  in  non-luminous  flame,  Eq.  (7-1)  may  be  inte- 
grated in  polar  spherical  coordinates  to  give 
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The  total  axial  component  of  radiative  flux  Is  then 
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(10-1) 


where  the  factor  coT  represents  the  black  body  radiation  at  temperature 
T,  and  terms  enclosed  In  brackets  are  the  effective  emlsslvlty  of  the 
fluid,  e(T,  pL) . Thus,  in  order  to  evaluate  Eq.  (10-1),  experimental 
measurements  of  e(T,pL)  as  a function  of  T and  pL  for  the  fire  environ- 
ment must  be  made. 

Another  radiative  transfer  problem  that  appears  to  be  applicable  in 
systems  analysis  of  fires  is  that  of  the  "stationary  radiation  front." 

If  conditions  permit,  the  temperature  vs.  x (the  horizontal  coordinate) 
can  be  related.  The  simple  geometry  of  an  infinite  flat  plane  with  a 
moving  fire  front  of  velocity  v^  at  tQ  and  a constant  specific  heat  cv 
is  considered.  The  kinetic  energy  and  thermal  conductivity  of  the  fluid 
are  neglected,  and  external  forces  are  assumed  to  be  zero.  The  absorption 
is  assumed  to  be  so  intense  that  the  radiation  flux  can  be  expressed  in 
terms  of  the  Rosseland  mean  free  path  as 


pv  c (T  - T ) ■ q„  - (qD  ) (11-1) 

xv  » Rx  Rx® 

It  is  assumed  that  T^  and  (q^)*,  are  much  smaller  than  the  correspond- 
ing values  at  the  fire  front.  For  the  above  assumptions,  the  energy  balance 
for  the  steady  state  can  be  represented  as 
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Differentiation  of  (11-1)  yields  — for  use  in  (12-1). 

ax 

The  mean  free  path  is  assumed  to  vary  as  the  s power  of  temperature; 
i.e.,  XR  - AqT8. 

A solution  of  Eq.  (12-1)  is  the  following  [see  reference  (4)]: 
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The  above  account  of  radiation  in  a fire  environment  iu  perhaps 
overly  simplified.  However,  a solution  to  the  integro-dlfferentlal 
equation  for  radiation  transfer  is  not  practical  at  this  time. 

THE  THERMAL  CONVECTION  MODEL 

For  the  work  reported  herein,  a model  for  convection  was  used  that 
Ignores  radiative  transfer  and  treats  a fluid  consisting  only  of  air, 
water  vapor,  and  condensed  drops.  In  this  model  the  dynamics  is  emphasized 
rather  than  the  chemistry.  Cloud  microphysics  is  treated  parametrically 
rather  than  explicitly. 

Thermal  convection  on  the  scale  of  a fire  situation  was  simulated 
by  numerical  solution  of  the  Boussinesq  equations  of  motion,  the  first 
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law  of  thermodynamics , and  the  equations  of  continuity  for  air  and 
water  substance.  The  Boussinesq  approximation  describes  a fluid 
that  is  incompressible  but  not  hydrostatic.  Hence,  a stream  function 
ip  can  be  defined  such  that  the  horizontal  and  vertical  components  of 
the  wind  are 
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in  two-dimensional  flow  in  Cartesian  (x,y,z)  coordinates,  and 
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in  axially  symmetric  flow  in  cylindrical  (r,0,z)  coordinates. 

A vorticity  equation  derived  from  the  equation  of  motion  is  used 
to  compute  the  changes  in  and  consequently  in  the  components  of  v. 
Concurrently,  the  first  law  of  thermodynamics,  the  equation  of  continuity 
for  water  substance,  and  a parameterization  of  the  mlcrophyslcs  of  con- 
densation and  evaporation  are  used  to  determine  changes  in  temperature 
departure  and  moisture  content;  hence  in  buoyancy. 

This  model  is  based  on  that  developed  by  Murray,  et  al.  (5,6,7)  to 
simulate  cumulus  convection  and  bears  some  similarity  to  the  cumulus 
model  of  Orville  (8).  The  boundary  conditions  require  that  the  normal 
component  of  wind  at  the  lower,  outer,  and  upper  boundaries  vanish. 

The  outer  and  upper  boundaries  are  chosen  at  a distance  from  the  disturbed 
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reglon  sufficient  to  minimize  unwanted  interactions.  In  addition, 

the  horizontal  component  of  wind  and  Lhe  horizontal  derivative  of 

the  vertical  component  of  wind  are  taken  to  be  zero  at  the  axis  of 

symmetry.  Eddy  diffusion  of  the  Flckian  type  is  included  in  the 

equations.  It  has  been  found  in  the  present  instance  that  eddy-diffusion 

2 -1 

coefficients  of  200,  40,  and  500  m sec  for  momentum,  moisture,  and 
temperature  respectively  lead  to  values  of  those  properties  that  are 
of  the  order  of  magnitude  of  observed  values. 

In  the  k isent  experiment  the  rectilinear  geometrv  was  used,  and 
the  convection  was  driven  by  a time-dependent  temperature  distribution 
specified  at  the  earth's  surface.  This  forcing  function  was  based  on  ob- 
servations of  a typical  forest  fire  (9);  it  is  illustrated  in  Fig.  1. 

The  computed  wind  fields  at  one,  two,  four,  and  six  minutes  of 
simulated  time  are  shown  in  Figs.  2-5.  The  arrows  fly  with  the  wind, 
and  each  half  barb  represents  a speed  of  5 m sec  , each  full  barb 
10  m sec  \ and  each  pennant  50  m sec  X.  The  direction  and  speed  .’..e 
also  given  by  the  digits  near  the  heads  of  the  arrows.  The  first  two 
digits  (or  one  digit  with  the  leading  zero  dropped)  show  the  direction 
in  tens  of  degrees  from  which  the  wind  is  blowing,  measured  clockwise 
from  the  "12  o'clock"  position.  The  last  two  numbers  give  the  speed  in 
m sec  The  values  shown  are  found  to  be  in  agreement  with  observations. 

The  departures  of  temperature  (#C)  from  the  initial  value  are  con- 
toured in  Figs.  6-10.  The  occasional  crossed  contours  were  Introduced 
by  Interpolation  made  by  the  automatic  contouring  program.  It  will  be 
noted  that  by  6 minutes,  the  excess  of  temperature  in  the  core  of  the 
fire  is  greater  than  100°C.  Unfortunately,  there  appears  to  be  little 
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experimental  data  available  on  temperature  departures  in  the  vicinity 
of  large  fires,  and  consequently  it  is  difficult  to  relate  the  simu- 
lated results  to  reliable  measurements. 

In  its  present  state  of  development,  the  program  can  handle  two 
geometries:  cylindrical  coordinates  with  axial  symmetry  and  Cartesian 

coordinates  with  symmetry  about  the  Y-Z  plane.  The  version  in  cylindrical 
coordinates  cannot  realistically  simulate  ambient  winds,  but  can  only 
simulate  an  Isolated  fire  under  calm  conditions.  It  is  capable  of  por- 
traying the  low-level  Indraft,  the  high-level  outdraft,  the  central  up- 
draft, and  the  surrounding  subsidence. 

The  version  in  Cartesian  coordinates  has  tie  disadvantage  of  pro- 
ducing downdrafts  that  are  somewhat  stronger  than  those  observed.  More- 
over, with  the  symmetry  condition  imposed,  it  can  handle  ambient  winds 
no  better  than  the  version  in  cylindrical  coordinates.  However,  it 
should  be  an  easy  matter  to  remove  the  symmetry  condition,  whereupon  the 
Cartesian  version  could  readily  take  care  of  an  ambient  wind  with  or 
without  vertical  shear,  not  to  mention  a non-level  lower  boundary.  Thus 
it  would  be  suitable  for  modeling  a firefront. 

Still  other  modifications  to  increase  the  utility  of  this  program 
for  simulating  the  dynamics  of  a real  fire  and  its  associated  environ- 
ment are  needed  and  feasible.  Besides  those  already  mentioned,  the  pro- 
vision for  modeling  fire  spread  by  winds  and  firebrands  and  the  provision 
for  modeling  the  spread  of  noxious  gases  come  readily  to  mind. 

The  model  is  sufficiently  general  and  flexible  that  all  of  these 
Improvements,  and  others  as  well,  can  readily  be  incorporated.  The  greatest 
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difficulty  foreseen  is  overcoming  the  numerical  instability  that  results 
from  the  strong  temperature  gradients  in  a fire.  It  appears,  however, 
that  with  a fine  mesh,  sufficiently  short  time  step,  and  careful  choice 
of  numerical  techniques,  this  problem  can  be  overcome.  With  the  model 
further  developed  as  described  above,  numerical  simulation  of  fires  and 
their  associated  environments  for  more  realistic  conditions  holds  promise 
of  becoming  a useful  analytic  tool. 
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FORCING  FUNCTION  EMPLOYED  TO  CAUSE  NUMERICAL  SIMIII  ATTON 
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Fig.  2— Wind  field  at  1 min  of  simulated  time. 
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Fig.  3— Wind  field  at  2 min  of  simulated  time. 
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Fig.  4--W1nd  field  at  4 min  of  simulated  time. 
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Fig.  5— Wind  field  at  6 min  of  simulated  time 


Fig.  6.— Temperature  contours  at  1 min. 
of  simulated  time. 

Scale:  1"  s 343  meters 


c%0 


Fig.  7.— Temperature  contours  at  2 min. 
simulated  time 
Scale:  1”  s 343  meters 
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